Source code for hysop.backend.device.codegen.structs.mesh_info

# Copyright (c) HySoP 2011-2024
#
# This file is part of HySoP software.
# See "https://particle_methods.gricad-pages.univ-grenoble-alpes.fr/hysop-doc/"
# for further info.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.


import numpy as np

from hysop.tools.htypes import check_instance, first_not_None
from hysop.tools.misc import upper_pow2_or_3
from hysop.backend.device.codegen.base.enum_codegen import EnumCodeGenerator
from hysop.backend.device.codegen.base.struct_codegen import StructCodeGenerator
from hysop.backend.device.opencl.opencl_types import OpenClTypeGen
from hysop.constants import BoundaryCondition

BoundaryConditionEnum = EnumCodeGenerator(
    BoundaryCondition, comments="Mesh boundary conditions enumeration"
)


[docs] class MeshBaseStruct(StructCodeGenerator): def __init__(self, typegen, vsize, typedef=True): assert vsize in [1, 2, 3, 4, 8, 16] name = f"{typegen.fbtype[0]}MeshBase{vsize}D" dtype, comments = MeshBaseStruct.build_dtype(typegen, vsize) if isinstance(typedef, bool): if typedef is True: typedef = name + "_s" else: typedef = None super().__init__( name=name, dtype=dtype, typegen=typegen, typedef=typedef, comments=comments ) self.reqs["boundary_enum"] = BoundaryConditionEnum self.vsize = vsize
[docs] @staticmethod def build_dtype(typegen, vsize): tg = typegen intn = tg.dtype_from_str(f"int{vsize}") floatn = tg.dtype_from_str(f"fbtype{vsize}") dtype = [] for fn in ["resolution", "compute_resolution"]: field = (fn, intn) dtype.append(field) for fn in ["xmin", "xmax", "size"]: field = (fn, floatn) dtype.append(field) dtype.append(("lboundary", intn)) dtype.append(("rboundary", intn)) comments = [ "Resolution of the mesh -including- ghosts", "Resolution of the mesh -excluding- ghosts", "Coordinates of the 'lowest' point including ghosts", "Coordinates of the 'highest' point including ghosts", "size = xmax - xmin", "Left boundary conditions", "Right boundary conditions", ] return np.dtype(dtype), comments
[docs] def create( self, name, resolution, compute_resolution, boundaries, xmin, xmax, size, **kargs, ): vsize = self.vsize dtype = self.dtype def extend(var, d=0): return np.asarray(tuple(var) + (d,) * (vsize - len(var))) tg = self.typegen lboundary, rboundary = boundaries lboundary = extend(lboundary, BoundaryCondition.NONE) rboundary = extend(rboundary, BoundaryCondition.NONE) _lboundary = [bd() for bd in lboundary] _rboundary = [bd() for bd in rboundary] mesh_base_vals = { "resolution": tg.make_intn(resolution, vsize), "compute_resolution": tg.make_intn(compute_resolution, vsize), "lboundary": tg.make_intn(_lboundary, vsize), "rboundary": tg.make_intn(_rboundary, vsize), "xmin": tg.make_floatn(xmin, vsize), "xmax": tg.make_floatn(xmax, vsize), "size": tg.make_floatn(size, vsize), } var = np.empty(shape=(1,), dtype=dtype) for k in mesh_base_vals: var[k] = mesh_base_vals[k] lboundary = BoundaryCondition.array_variable( "lboundary", typegen=tg, vals=lboundary ) rboundary = BoundaryCondition.array_variable( "rboundary", typegen=tg, vals=rboundary ) value = dict( resolution=extend(resolution), compute_resolution=extend(compute_resolution), xmin=extend(xmin), xmax=extend(xmax), size=extend(size), ) var_overrides = {"lboundary": rboundary, "rboundary": lboundary} cg_var = self.build_codegen_variable( name=name, value=value, var_overrides=var_overrides, **kargs ) return (var, cg_var)
[docs] class MeshInfoStruct(StructCodeGenerator): def __init__(self, typegen, vsize, typedef=True, mbs_typedef=None): assert vsize in [1, 2, 3, 4, 8, 16] name = f"{typegen.fbtype[0]}MeshInfo{vsize}D" bname = f"{typegen.fbtype[0]}MeshBase{vsize}D_s" mbs_typedef = first_not_None(mbs_typedef, bname) if isinstance(typedef, bool): if typedef is True: typedef = name + "_s" else: typedef = None dtype, comments, ctype_overrides, reqs = MeshInfoStruct.build_dtype( typegen, vsize, mbs_typedef=mbs_typedef ) super().__init__( name=name, dtype=dtype, typegen=typegen, typedef=typedef, comments=comments, ctype_overrides=ctype_overrides, ) for req in reqs: self.require(req.name, req) self.mesh_base = reqs[0] self.vsize = vsize
[docs] @staticmethod def build_dtype(typegen, vsize, mbs_typedef=None): tg = typegen mbs_typedef = first_not_None( mbs_typedef, f"{typegen.fbtype[0]}MeshBase{vsize}D_s" ) intv = tg.dtype_from_str(f"int{vsize}") floatv = tg.dtype_from_str(f"fbtype{vsize}") i = 0 dtypes = [] def _append(dtype): dtypes.append(dtype) _append(("dim", np.int32)) i += 1 for fn in ["start", "stop", "ghosts"]: _append((fn, intv)) i += 1 for fn in ["dx", "inv_dx"]: _append((fn, floatv)) i += 1 mesh_base = MeshBaseStruct(tg, vsize, typedef=mbs_typedef) _append(("local_mesh", mesh_base.dtype)) _append(("global_mesh", mesh_base.dtype)) i += 2 comments = [ "Dimension of the mesh", "Index of the first local compute point in the global grid", "Index of the last local compute point in the global grid", "Number of ghosts in each direction", "Space discretization", "1/dx", "Local mesh", "Global mesh", ] return dtypes, comments, None, [mesh_base]
[docs] def create( self, name, dim, start, stop, ghosts, dx, local_mesh, global_mesh, **kargs ): vsize = self.vsize if dim > vsize: msg = f"Dim should be less or equal to {vsize}, got dim={dim}." raise ValueError(msg) tg = self.typegen vsize = self.vsize dtype = self.dtype dx = np.asarray(dx) mesh_info_vals = { "dim": tg.make_intn(dim, 1), "start": tg.make_intn(start, vsize), "stop": tg.make_intn(stop, vsize), "ghosts": tg.make_intn(ghosts, vsize), "dx": tg.make_floatn(dx, vsize), "inv_dx": tg.make_floatn(1.0 / dx, vsize), "local_mesh": local_mesh[0], "global_mesh": global_mesh[0], } def extend(var, d=0): if isinstance(var, np.ndarray): _var = np.full(shape=(vsize,), fill_value=d, dtype=var.dtype) _var[: var.size] = var return _var else: return np.asarray(tuple(var) + (d,) * (vsize - len(var))) value = dict( dim=dim, start=extend(start), stop=extend(stop), ghosts=extend(ghosts), dx=extend(dx), inv_dx=extend(1.0 / dx), ) local_mesh[1].name = "local_mesh" global_mesh[1].name = "global_mesh" var_overrides = dict( local_mesh=local_mesh[1], global_mesh=global_mesh[1], ) cg_var = self.build_codegen_variable( name=name, value=value, var_overrides=var_overrides, **kargs ) var = np.empty(1, dtype) for k in mesh_info_vals: var[k] = mesh_info_vals[k] return (var, cg_var)
[docs] @staticmethod def create_from_mesh(name, typegen, mesh, **kargs): from hysop.mesh.cartesian_mesh import CartesianMeshView check_instance(name, str) check_instance(mesh, CartesianMeshView) check_instance(typegen, OpenClTypeGen) tg = typegen dim = mesh.dim vsize = upper_pow2_or_3(dim) btd = f"{tg.fbtype[0]}MeshBase{vsize}D_s" itd = f"{tg.fbtype[0]}MeshInfo{vsize}D_s" mesh_base = MeshBaseStruct(tg, vsize, typedef=btd) mesh_info = MeshInfoStruct(tg, vsize, typedef=itd, mbs_typedef=btd) start = mesh.global_start[::-1] stop = mesh.global_stop[::-1] dx = mesh.space_step[::-1] ghosts = mesh.ghosts[::-1] gresolution = mesh.grid_resolution[::-1] gcompute_resolution = mesh.grid_resolution[::-1] gxmin = mesh.global_origin[::-1] gxmax = mesh.global_end[::-1] gsize = mesh.global_length[::-1] gboundaries = mesh.global_boundaries[::-1] lresolution = mesh.local_resolution[::-1] lcompute_resolution = mesh.compute_resolution[::-1] lxmin = mesh.local_origin[::-1] lxmax = mesh.local_end[::-1] lsize = mesh.local_length[::-1] lboundaries = (mesh.local_boundaries[0][::-1], mesh.local_boundaries[1][::-1]) const = kargs["const"] if ("const" in kargs) else False gmesh = mesh_base.create( name + "_global_mesh", gresolution, gcompute_resolution, gboundaries, gxmin, gxmax, gsize, const=const, ) lmesh = mesh_base.create( name + "_local_mesh", lresolution, lcompute_resolution, lboundaries, lxmin, lxmax, lsize, const=const, ) (var, cg_var) = mesh_info.create( name, dim, start, stop, ghosts, dx, lmesh, gmesh, **kargs ) return (var, cg_var)
[docs] def build_codegen_variable(self, name, var_overrides=None, **kargs): tg = self.typegen if var_overrides is None: const = kargs["const"] if ("const" in kargs) else False var_overrides = dict( local_mesh=self.mesh_base.build_codegen_variable( "local_mesh", const=const ), global_mesh=self.mesh_base.build_codegen_variable( "global_mesh", const=const ), ) return super().build_codegen_variable( name=name, var_overrides=var_overrides, **kargs )
if __name__ == "__main__": from hysop.backend.device.codegen.base.opencl_codegen import OpenClCodeGenerator from hysop.backend.device.codegen.base.test import _test_typegen tg = _test_typegen("double", float_dump_mode="dec") vsize = 4 mbs = MeshBaseStruct(tg, typedef=f"{tg.fbtype[0]}MeshBase{vsize}D_s", vsize=vsize) mis = MeshInfoStruct( tg, typedef=f"{tg.fbtype[0]}MeshInfo{vsize}D_s", mbs_typedef=mbs.typedef, vsize=vsize, ) cg = OpenClCodeGenerator("test_generator", tg) # declare mesh MeshInfoStruct and its dependancies (MeshBaseStruct,MeshDirectionEnum,TranspositionStateEnum) cg.require("mis", mis) # create a local numpy and a codegen MeshInfoStruct variable local_mesh = mbs.create( "local", ( 10, 10, 10, ), ( 8, 8, 8, ), ((BoundaryCondition.PERIODIC,) * 3,) * 2, ( 0, 0, 0, ), ( 1, 1, 1, ), ( 1, 1, 1, ), ) global_mesh = mbs.create( "global", ( 100, 100, 100, ), ( 80, 80, 80, ), ((BoundaryCondition.PERIODIC,) * 3,) * 2, ( 0, 0, 0, ), ( 10, 10, 10, ), ( 10, 10, 10, ), ) (np_mis, cg_mis) = mis.create( "mesh_info", 3, ( 0, 0, 0, ), ( 1024, 1024, 1024, ), ( 1, 1, 1, ), (0.1, 0.2, 0.3), local_mesh, global_mesh, storage="__constant", ) # declare and intialize the nested struct in the __constant address space cg.jumpline() cg_mis.declare(cg) # show generated code in your editor cg.edit() # build generated code on every OpenCL device to check cg.test_compile()